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Abstract: We discuss the contribution of OZI-suppressed diagrams to the hyperfine split- 
ting of charmonium in lattice QCD. We study valence quark mass regions from strange to 
charm quark masses. No contribution of the disconnected diagram is seen in the vector me- 
son channel. In the pseudo-scalar channel and for valence quark masses around the strange 
quark, the disconnected contribution induces a considerable increase of the meson mass. 
This contribution quickly decreases as the quark mass increases. For charmonium the effect 
is very small although a decrease of the pseudoscalar mass induced by the disconnected 
contribution cannot be ruled out. 



Keywords: 



Lattice Quantum Field Theory, Lattice Gauge Field Theories, Lattice QCD 



Contents 



[l]. Introduction 



[2]. Formulation 

!T] Quarkonium correlation function 
Valence quark action 

|3|. Numerical simulations 
|3,1| Lattice setup 

Tuning of valence heavy quark action 
Connected correlators 
Evaluation of disconnected diagram 
Quarkonium correlators 



3.2 



3.3 








T5| 



EL Conclusions and discussion 



12 



1. Introduction 

Lattice QCD has been able to provide a wealth of results for the hadron spectrum, both in 
the quenched approximation and with dynamical quarks [Q]. Most of the low-lying hadron 
masses are reproduced under extrapolation to the continuum limit and to the physical u 
and d quark masses. Even in the quenched approximation there is a remarkable agreement 
with experimental data for the light hadron spectrum, deviations typically amounting 
to 10 per cent. The charmonium hyperfine splitting is, however, an exception to such 
success. Recent systematic computations in the quenched approximation have pointed out 
a significant discrepancy between the lattice results and experimental data: the former is 
30-40% smaller than the latter, AM = Mj^ - M Vc = 117 MeV. 

It is difficult to address the calculation of charmonium spectrum with current com- 
putational resources. Since the charm quark mass is not well below the accessible lattice 
cutoffs, lattice formulations that reduce the 0(am c ) lattice artifacts seem mandatory. The 
most promising, while brute force approach, is to use a relativistic formulation with suffi- 
ciently small lattice spacing and 0(a)-improved quark actions. This is the approach taken 
in Ref. Q by the QCD-TARO Collaboration, using the nonperturbatively 0(a)-improved 
Wilson quark action on quenched isotropic lattices Q. Numerical simulations with lattice 
cutoffs ranging from 2 to 5 GeV found AM = 77(2)(6)MeV in the continuum limit. 

Other approaches involve the use of effective heavy quark actions. Among them, non- 
relativistic QCD (NRQCD) has been investigated most extensively. The latest NRQCD 
quenched result is AM = 55(5) MeV (with the scale set by the P-S splitting) |3]. In this 
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Figure 1: Connected (Left) and OZI-suppressed (Right) diagrams contributing to the pseudoscalar 
(r = 75) and vector (r = 7 M ) channels. 



case, lattice artifacts are difficult to control since brute force elimination by taking the 
continuum limit is impossible. Relativistic formulations, such as the Fermilab approach 
H|, |H| and anisotropic lattices || |7[ ||, [J, have an advantage in this sense. However, 
they give results for the, continuum extrapolated, hyperfine splitting which also definitely 
deviate from experiment. Here it should be noted that for a quark mass not sufficiently 
smaller than the spatial lattice cutoff, the dynamics inside heavy quarkonia may not be 
precisely described by these actions whose spatial derivative terms are defined only with 
nearest neighbouring sites. This is because of the 0((a a p) 2 ) errors, where a a is the spatial 
lattice spacing and p the typical quark momentum. This error might become important 



for heavy quarkonium, for which p ~ am, J9J, 10 1 . This is in contrast with the situation for 



light and heavy-light hadrons, where p ~ ^-qcd- Therefore advantages of these relativistic 
formulations applied to heavy quarkonia are rather limited. 

From all the previous studies, disagreement between the quenched lattice calculation 
and experiment has been established for the charmonium hyperfine splitting. 

There are two candidates to explain this discrepancy. The first one is dynamical quark 
effects. Although a systematic study of these effects, involving a continuum extrapolation, 
has not been performed yet, several groups have tried to estimate them [11, 13]. A recent 



lattice computation with 2+1 flavours of improved staggered quarks at a -1 ~ 1.6 GeV has 
reported a hyperfine splitting AM = 97(2) GeV [JDJ. This value is still 20% smaller than 
the experimental one. Since they applied the Fermilab action with tadpole-improved tree- 
level value for the clover coefficient, the remaining discrepancy may be attributed to the 
O(cra) and 0((ap) 2 ) systematic errors. Systematic studies with higher lattice cutoff's and 
involving continuum limit extrapolations are strongly desired. 

Another possible contribution, which has not been incorporated in any of the lattice 
computations (quenched or unquenched) performed up to now, comes from OZI-suppressed 
(disconnected) diagrams as those in Fig. [j] |2|. Such diagrams must be included in the 
evaluation of correlators of unflavoured mesons, such as r\ c and J/ip. Although, according 
to the perturbative picture, this contribution is expected to be small in heavy quarkonium, 
it might be non-negligible compared to the, also small, hyperfine splitting. Moreover, its 
effect might be enhanced. This happens, indeed, in the light quark mass region, where the 
contribution of the disconnected diagram to the pseudoscalar channel is strongly enhanced 
by the ?7a(1) anomaly . It is therefore important to quantify the size of such contribution 



1 Note, however, that the contribution of the Ua(1) anomaly raises the mass of the pseudoscalar while 
not affecting that of the vector. This effect would induce a decrease of the hyperfine splitting, instead of 
the increase required to match the experimental value for charmonium. 
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to the charmonium correlator. 

The goal of this paper is to examine the second possibility. A similar analysis has 
recently been performed by McNeile and Michael (UKQCD Collaboration) in Ref. p0[ . 
We aim here at an exploratory study to estimate whether the disconnected diagram can 
give a significant contribution to the charmonium hyperfine splitting. This issue is not 
completely disentangled from the one of dynamical quark effects. Indeed, these effects can 
be particularly large for disconnected diagrams. In the quenched approximation closed 
quark loops, as those in the right panel of Fig. |], can only be connected through gluonic 
contributions. For unquenched simulations, however, they can also be connected by inser- 
tion of virtual quark loops. This induces a very different asymptotic behaviour between 
quenched and unquenched disconnected correlators, as discussed at length for the case of 
the lattice determination of the rf mass [|lj]]. In addition, through the disconnected di- 
agram charmonium states mix with glueballs and lighter quarkonium. Mixing with the 
latter is completely neglected in the quenched calculation, while mixing with the former is 
only partially taken into account. 

We do not intend here to analyze the failure of the quenched approximation in the 
calculation of the hyperfine splitting but to address the more general question of whether 
disconnected correlators can give a sizable contribution to this splitting. Since it is ex- 
pected that such contribution will be quickly obscured by statistical noise, large statistics 
is essential for the present study and this has dictated our choice of lattice parameters 
and lattice action. The number of configurations collected for the quenched calculation of 
the hyperfine splitting in || is not sufficient for this study and increasing the statistics for 
these large lattices is beyond our computational resources. For this reason we use, instead, 
a large set (3200) of available configurations with a rather coarse 12 3 x 24 lattice, with 
fixed lattice cutoff a -1 ~ 1.2 GeV and with two dynamical flavours of staggered quarks 
with am seSu = 0.10 (no attempt at a continuum extrapolation will be presented in this 
paper). Given the rather large sea quark masses adopted, we expect the effect of dynam- 
ical quarks to be relatively small, hence our results give also an insight for the quenched 
situation. It would be important to study the sea quark mass dependence to examine how 
such contribution shifts the charmonium masses for physical sea quark masses and Nf = 3. 

In order to treat charm quarks on our rather coarse lattice, we adopt the Fermilab 
quark action. Since the disconnected diagram contribution to the charmonium hyperfine 
splitting is hard to detect because of large statistical noise, we also compute this contri- 
bution for lighter valence quarks, and estimate how it varies as the valence quark mass 
is increased towards that of the charm quark. This helps us determine its value at the 
charm quark mass. This is the only purpose of our varying the valence quark mass. For 
lighter valence quarks, the inconsistency between our sea and valence quark actions would 
generate systematic errors, which would cause additional problems for the determination 
of the meson spectrum and of the hyperfine splitting. We do not consider this regime here. 

The paper is organised as follows. In the next section we describe the set up for our 
calculation. Section [3| presents our results. The last section is devoted to conclusions and 
discussion. 
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2. Formulation 



2.1 Quarkonium correlation function 

We define the quarkonia correlation function with a quarkonium operator 0(x,t), 

0{x,t) =Y,q{x + y,t)Tq(x,t)v(y), (2.1) 
y 

where <p(y) is a smearing function and T is 4 x 4 matrix which specifies the quantum 
numbers of the quarkonium state. 

We take T = 75 and 7^ for the pseudoscalar and vector channels, respectively. Since 
the operator 0(x,t) contains quark fields of the same flavour, by contracting the quark 
lines the correlator decomposes into a connected C con (t) and a disconnected part Cdi S (t), 

Cfuiiit) = ^(0^)0(0,0)) = C con {t)+C dis {t). (2.2) 

X 

Using the quark propagator D~ l (x, t; x', t'), where D is the Dirac operator, we can write 
C con (t) and C dis (t) as follows, 

C C on{t) = Y1 (Tr^(y)D- lt (x + y,t;0,0)7 5 r^(z)D- 1 (x,t;i',0)r75]) , (2.3) 

x,y,z 

C dls (t) = i- (L(t)*L(0)) , (2.4) 

where TV is the trace over the colour and spinor indices, V3 the spatial volume. The quark 
loop diagram L(t) is defined as 

L(t) = Y,Tr[i P (y)D- 1 (x,t;x + gt)T] . (2.5) 

%,y 

As the valence quark mass increases, the disconnected correlator is expected to decrease 
compared to the connected part and the extraction of a signal for the disconnected part may 
become increasingly difficult. It is therefore crucial to improve the quarkonium operator 
so as to increase the overlap with the ground state. For this purpose, we use spatially 
extended operators with smearing function (p(x) in the Coulomb gauge. For (p{x) we adopt 
Gaussian functions with width around the expected charmonium radius and select the best 
one among them. 



In order to evaluate the trace in Eq. ( |2.5D , we employ the complex Z2 noise method [16]. 
This method is a popular technique to evaluate quark loop contributions to correlators. 
An application to the smeared operator, Eq. ( |2.5| ), is straightforward. We note that one 
needs to solve the quark propagator only once for each noise vector for all the smearing 
functions applied. 

2.2 Valence quark action 

The coarse lattice we employ does not allow for an isotropic formulation of heavy valence 
quarks whose mass is not sufficiently smaller than the lattice cutoff. For this reason we 
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adopt the Fermilab action 



xy 



(l-74)U^x)6 x+iy + (l+^)U 4 (x - 1)5. 



x—A,y 



x.y 



-«<r£i [0 - li)Ui{x)5 x+ly + (r + 7i)0i(s-*)£ B _i, y j 

(x,y)S x , y }q(y), 



(2.6) 



where the spatial Wilson parameter is set to r = 1. The parameter 7^ is to be tuned so 
that the rest mass, M\ = E(p = 0), equals the kinetic mass, 



1 

M~ 2 



d 2 E(p) 



dp 2 



(2.7) 



p=0 



for, for example, a meson dispersion relation. We define k by incorporating tadpole im- 
provement [|l7]] as 

- = — -2{ lF + 3r-4) (=2(m + 4)), (2.8) 

K Uq K a 

where niQ is the bare quark mass [jOl]. As the mean-field value of link variable, no, we 
adopt the average value in the Landau gauge, no = (Trf7 M (x)}/3. The value of jf is tuned 
for each value of k. Then the chiral extrapolation is to be performed in 1/k, a step that 
is not taken in this paper. The clover coefficients, ce and ce, control the 0{a) systematic 
uncertainty. In this work, we adopt the tadpole improved values [17], 

c E = c B = 



(2.9) 



This action allows, in principle, for a relativistic treatment of heavy quarks but, as 
pointed out in Sec. |l], even with nonperturbatively tuned jf, the results suffer from 
0((ap) 2 ) errors, where p ~ am q in the case of heavy quarkonium, in addition to the 
O(aa) error from the clover terms. We expect this effect not to be essential for a qual- 
itative estimate of the size of the contribution from the disconnected diagram. Better 
control over these errors is required, however, for a quantitative evaluation and continuum 
extrapolation of the contribution. 



3. Numerical simulations 
3.1 Lattice setup 

The numerical simulation is performed on lattices of size 12 3 x 24, with two flavours of 
staggered dynamical quarks. The gauge action is the standard Wilson plaquette action 
with = 5.50. The dynamical staggered quark mass is m sea a = 0.10. The configuration 
is generated with the hybrid R- algorithm with 5t = 0.02 and unit length of trajectory. We 
prepare 32 independent initial configurations and generate configurations in parallel. Each 
configuration is separated by 5 trajectories, after 900 trajectories for thermalization. 

The lattice cutoff scale is set by Sommer's hadronic radius, rp, which is defined through 
r^Firo) = 1.65, where F(r) is the force between static quark and antiquark [19]. Setting 
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Figure 2: Tuned parameter 7^ using the meson dispersion relation. 
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Table 1: Parameters of the valence quarks. N con f is the number of gauge configurations used; 
they are generated from 32 independent initial configurations. The pseudoscalar and vector meson 
masses are extracted from connected correlators with smeared operators. 



1.1 

1 

> 0.9 
0.8 
0.7 



ro = 0.5 fm yields a -1 = 1.2012(60) GeV. For tadpole improvement [17], the mean-field 
value is determined in the Landau gauge as uq = (TrC/ /i (x))/3, giving no = 0.824285(89) 
and clover coefficient csw = 1-785 for our lattice. In the following analysis, the statistical 
error is estimated by the jackknife method. 



3.2 Tuning of valence heavy quark action 

In this subsection, we describe the tuning of the parameter 7^ in the valence quark action. 
The tuning procedure is the same as for anisotropic lattices [18] . At each hopping parameter 
k, the coefficient 7^ in the quark action (|2.6| ) is determined nonperturbatively using the 
meson dispersion relation. For the meson dispersion relation, we assume the relativistic 
form 

E(p)* = M? + ^-p 2 + 0(p% (3.1) 
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where Mi and M2 are the rest and kinetic masses of the meson. The value of 7^ is tuned 
so that Mi = M 2 holds. 

The tuning of jf is performed with connected correlators of point operators, namely 
<p(x) = 5(x), on 400 configurations. The meson energies are fitted to a quadratic form 
for the pseudoscalar and vector channels. The obtained £p = y 1 'M2/M1, the fermionic 
anisotropy, are spin averaged and fitted to a linear function of jf- Interpolating to 
(M1/M2) = 1, the tuned value of 7_f, 7^, is determined. 

The result for 7^ is displayed in Fig. |2[ Since in the light quark mass region the 
Fermilab action smoothly tends to the standard clover quark action, the value of 7^ should 
approach unity as the quark mass decreases. From a tree level analysis, it should be a 
decreasing function of 1/k, and in the light quark mass region, a linear dependence in 
is expected (m q = (1/k — 1/k c )/2, where k c is the critical hopping parameter (which is 
not determined on the present lattice). In Fig. although the decreasing tendency as 1/k 
increases is indeed observed, the quark mass dependence in the small mass region seems not 
in agreement with the expected behaviour. This is presumably due to large lattice spacing 
artifacts, one of which apparently comes from the assumed form of the meson dispersion 
relation. This uncertainty may also cause 7J, to approach a value slightly different from 
unity J18| |. Although we will not try to correct these systematic errors for the present 
qualitative estimate of the hyperfine splitting, it is essential to have them under control in 
order to extract reliable quantitative determinations. 

The value of k corresponding to the charm quark mass is determined by interpolating 
the results for the vector meson mass so as to reproduce the physical J/ip meson mass. 
The value of 7^ is also interpolated to that k value. The resulting k and 7^ are listed in 
the last line of Table |], and displayed as the rightmost point in Fig. ||| together with the 
other five cases listed in Table |l[ 

3.3 Connected correlators 

Charmonium correlators are computed for the quark parameters listed in Table |l[ Let us 
start with the connected correlator, Eq. (|2.3| ). In the following analysis, the error due to 
the uncertainty in 7^ is not evaluated. 

We first observe the efficiency of the smearing technique on the connected correla- 
tors. Figure || displays the effective mass plot for the connected correlators with local and 
smeared operators in pseudoscalar and vector channels. Here the effective mass is defined 
with 

C(t) cosh {{T/2-t)m eff {t)} 

C(t + 1) cosh[(T/2-t-l)m e //(i)]' 1 ' ' 

We use smearing functions with Gaussian form and several width values, and select the one 
that gives better plateaus in the effective mass plot. Figure || clearly shows that the smeared 
operator considerably enhances the overlap with the ground state, which dominates the 
connected correlator beyond t = 4. 

In Table |l| we list the meson masses extracted from the connected correlators. At the 
charm quark mass, the obtained hyperfine splitting, AM ~ 81 MeV, is consistent with 
previous results 13], and again far below the experimental value. 
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Figure 3: Effective masses of connected correlators with local and smeared operators. 
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Figure 4: Typical samples of the dependence on the number of noise vector of pseudoscalar 
disconnected correlators at t — 3. Since the results are averaged over configurations, the error 
includes both that of the noise method and that of the ensemble average. 



3.4 Evaluation of disconnected diagram 



For an evaluation of the disconnected correlators, Eq. (2.4), we apply the Z2 noise method. 
Figure || shows typical samples of the pseudoscalar disconnected correlator at a time slice 
t = 3 as a function of the number of noise vectors. The top and bottom panels display 
the results for k = 0.11013 and 0.10476. Since the results are ensemble averages over 
configurations, the displayed error includes both the quantum fluctuations and the error 
due to a finite number of noise vectors. 

In order to estimate the contribution of disconnected diagrams we have to control 
both errors. A reasonable approach is to choose the number of noise vectors, N]yy, so as 
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to keep both errors at similar levels. Although we can find an "optimal Nnv" at each t 
by observing Nnv dependence of the error displayed in Fig. f|, it is sufficient to roughly 
estimate such Njyy, since increasing the number of configurations, N con f, also reduces the 
error of noise method and we choose rather large values of N con f. The number of noise 
vectors, Njyy, and gauge configurations, N con j selected for each valence quark mass are 
displayed in Table |l[ As the quark mass increases, the signal to noise ratio becomes worse 
rapidly. Thus a larger number of configurations as well as of noise vectors are prepared in 
this calculation. 

3.5 Quarkonium correlators 

In this section we evaluate the relevance of the disconnected part of the meson correlator 
for the meson spectrum. Figs. |B| and |6| show the ratio R(t) = Cdi s {t)/C con {t). R(t) is often 
used to evaluate the contribution of the disconnected diagram [jiq| . If the ground state 
dominates both the connected and full correlators, their asymptotic time dependence can 
be described by single exponentials: 

Ccon(t) = A con exp (—m con t) , (3-3) 
Cfuii{t) = C con {t) + C dis (t) = Af u u exp (-mf u ut). (3.4) 

Then the ratio behaves as 

R (t) = 77^4^ = ~~T^~ ex P {-( m full ~ m con )t} - 1. (3.5) 

^con\t) -^-con 

This ratio is useful for exploring the sign and magnitude of the mass difference {jnt u u — 
m con) when the signal is noisy. In Figs. [5] and ^, we plot both point and smeared correlator 
ratios R(t). For the latter, effective masses for full and connected correlators are also 
displayed in the right panels of Figs. |5| and ||. 

For the vector channel, we find no sizable contribution from the disconnected diagram 
in the whole quark mass region explored. Similar results have been also reported in ||15|| . 
This is consistent with the OZI suppression. In contrast, in the pseudoscalar channel, a 
clear signal is observed in the light quark mass region giving rise to an increase of the 
pseudoscalar mass. This is consistent with previous works on the flavour singlet r/ (r/) 
meson and with theoretical expectations. It is important to note that this effect is opposite 
to what would be needed to match the experimental result for charmonium. Given that the 
vector mass receives no sizable contribution, a decrease of the r\ c mass by about 30-40 MeV 
is required to bring the lattice measurement of the hyperfine splitting up to the value of 
117 MeV measured in experiment. Indeed, the mass difference between full and connected 
correlators rapidly decreases as the quark mass increases, leaving in principle room for a 
change in the sign of the mass difference. In the left bottom panel of Fig. |(| the dashed 



curves represent R(t) from Eq. 3.5 for mf u u — m con = ±20 MeV, assuming Af u u/A con = 1. 
In the scale of the plot, a small deviation of Af u n/A con from one would induce a shift of 
the curves, keeping the slopes almost constant. Due to the large statistical errors for t > 2 
we cannot make any definite statement about the slope of the data: it is compatible with 
zero within errors but small slopes of order ±20 MeV cannot be ruled out. 
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Figure 5: Ratio of disconnected to connected correlators (left panels) and effective mass of full and 
connected correlators (right panels) for k = 0.11294 (top), 0.11013 (middle), and 0.10732 (bottom). 



To better settle this point, we plot in Fig. [7] the mass difference between the full and 
connected correlators for the pseudoscalar channel as a function of the vector meson mass. 
Here the mass difference is defined by the difference between effective masses of full and 
connected correlators at t = 1 and t = 2. At these time slices, the effective mass plots 
in Figs. H and ^ do not exhibit clear plateau behaviour. However, the contribution of the 
excited states partially cancel in the difference. In fact, the results for t = 1 and 2 differ 
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curves represents the cases where m f u u — m con — ±20 MeV with A f u u /A con = 1 . 



from each other only slightly. We do not plot the result for t = 2 at the charm quark mass 
because the result is too noisy. Figure again shows that the mass difference between full 
and connected correlator is positive for light quark masses. It rapidly decreases as the quark 
mass increases, becoming almost zero at about half the charm quark mass. Our results at 
the charm quark mass are not conclusive. At the one sigma level they are compatible both 
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Figure 7: The mass difference between the fuii and connected correlators, in the pseudoscalar 
channel, defined with the effective masses at t = 1 and t = 2 as a function of vector meson mass. 
Wc call the effective mass at "t = 1" that derived from comparing t = 2 and t = 1. 



with zero mass difference or with a negative mass difference of about —20 MeV. This is 



consistent with the results found by McNeile and Michael (UKQCD Collaboration) in [20]. 



4. Conclusions and discussion 

We have investigated the contribution of disconnected diagrams to the hyperfine splitting 
of quarkonium in a quark mass region ranging from strange to charm quark masses. As 
displayed in Figures pi ||, we found almost no contribution of disconnected diagrams to 
the correlators in the vector channel in the whole quark mass region explored in this 
work. For the pseudoscalar channel, however, there is a sizable contribution around the 
strange quark mass region which quickly decreases as the quark mass increases and almost 
vanishes around half the charm quark mass. At the charm quark mass, the contribution of 
the disconnected correlator is very small, in agreement with the expected OZI suppression. 
Given our large statistical errors we cannot, however, rule out that the mass difference 
between full and connected correlators becomes negative, although small. In this respect 



we agree with the results previously found by McNeile and Michael [20] who found room 
for a mass difference of the order of —20 MeV . 

To determine how large is the contribution of disconnected diagrams for more physical 
situations, one needs to perform the chiral extrapolation in the sea quarks and the con- 
tinuum limit. As argued in the introduction, going to lighter sea quarks may considerably 
modify the contribution from disconnected diagrams. In particular, an effect that could 
induce a further decrease in the rj c mass would be mixing with a pseudoscalar glueball were 
this to be lighter than the r\ c . This is a possibility that has been discussed also in p0| . 
although it has been considered not very likely given the fact that lattice determinations 
of the pseudoscalar glueball gives masses below that of the r] c meson. 
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